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' Abstract. Fractional Brownian motion belongs to a class of long 

{/y . memory Gaussian processes that can be represented as linear func- 

' tionals of an infinite dimensional Markov process. This leads nat- 

£S| ■ uraly to : 

• An efficient algorithm to approximate the process. 

• An ergodic theorem which applies to functionals of the type 



CM 



x 



4>(Vh(s)) ds where V h (s) 



[ h{t-u)dB u 
Jo 



1. Introduction 



Fractional Brownian motion of Hurst parameter H e (0, 1) is a real 
centered Gaussian process {W^ 1 ; t > ) of covariance 

(1) E [W t H W s H ] = ^El(s 2H + t 2H -\t- s\ 2H ) (s, t > 0) . 

00 ■ 

. One easily sees that it starts from zero, Wq = almost surely, it has 

stationary increments (W£ s — ; t > ) = (W t H ; t > ) since 

(2) E[(W t H -W s H f] =c(H)\t-s\ 2H (s,t>0). 

Fractional Brownian motion enjoys the following scaling (self-similarity) 
property 



(3) (±W% ;t>0)±(W t H ;t>0) (a>0). 

(X 

The increments are independent if and only if H = | (Brownian mo- 
tion). If H < |, they are negatively correlated, and if H > | they are 
positively correlated. Furthermore it is a long memory process, since, 
for H ^ 2> the covariance between far apart increments decrease to 
zero as a power law: 

Cov{W H {t) - W H (0), W H (t + u)- W H (u)) 
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Fractional Brownian motion has been used as a model of signal or 



noise in various domains: geophysical data [10, p|, communication pro- 



cesses fL2| , |I1|| , see also the references therein. Therefore it is important 
to provide good numerical approximating schemes. 

Problem I: Approximation of fractional Brownian motion. We 

shall not linger on the first (naive) idea that may cross our mind. We 
simulate an n-dimensional Gaussian random vector with the same dis- 
ribution as (W£*, . . . , W^ n ) and we perform a linear or spline inter- 
polation. This amounts to computing the square root of the n x n 
covariance matrix and is computationnaly inneficient. 
A second solution is to use the representation introduced by Mandel- 
brot and Van Ness JL3|: if a = H + |, 

1 /"* 

rH 1 / IU „,\o-l / „.\a-l\ 



W» = — ((t - u)^ 1 - (-u)r ') dB u 
1 (a J J-oo 

where B is a standard Brownian motion on the real line. Without any 
loss in generality, we can say that the problem is to approximate the 
finite, but long, memory part: 

t 



V(t) = V h (t) = [ h(t-s)dB, 
Jo 



where (B t ; t > ) is a standard Brownian motion. The naive way to 
do this is to fix a small step |A n | = - for the partitions we use 

1 2 
A n = 0<S! = -<s 2 = -<•••■ 
n n 

Hence, if t = r/n, then V(t) is approximated by 

* - *) tow - a*) 

Si) 6 



fn ^ 



-L(h{~)S r + /i(2/n)e r -i + 
>n n 



where ^ = \/n(B s . +1 — B s .) are independent identically distributed 

standard normal random variables. 

For the next step, t + At = t + - = (r + l)/n, we have 

V(t + At) ~ -L(h(l/n)£ r+1 + h(2/n)£ r + ■ ■ •) 
<n 



Therefore, all we have to do is to compute the h(i/n), generate the 
and keep them all in memory for the next step. 

Unfortunately, if we want to keep the approximation error low (that is 
n big), this amounts to a lot of data to keep in memory. 
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The natural idea is to impose a threshold, that is to compute h(i/n) 
and generate £ only for 1 < i < M n . This may work well for smooth 
functions h, but not for 

T(a) 

where | < a < 1, since h explodes near zero, and thus our partition is 
not dense enough near t ! 

Unfortunately this example is important since it is in essence fractional 
Brownian motion of index H — a — |. 

In j7|,Carmona, Coutin and Montseny derived an approximation scheme 
from a representation formula inspired by the diffusive input-output 
model of fractional integration developped by Audounet, Montseny adn 
Mbodje ii. 

More precisely Fubini's stochastic theorem implies that if h has the 
"spectral" representation 

Mu) =/e-» W („>0) 

where \i is a positive measure on (0, +oo), finite on compact sets, and 
such that J(l A x~?) fi(dx) < +oo, then 

POO 

V h (t)= / Y t *n(dx) 
Jo 

where 

Y t x = ! e - x(t - s) dB s . 
Jo 

Observe that for fixed x > 0, (Yf ; t > ) is an Ornstein-Uhlenbeck 
process of parameter x, that is a Gaussian semi-martingale Markov 
process solution of the stochastic differential equation 

dY t x = dB t -x Y t x dt . 

Therefore, Vh = is a linear functional of the infinite (here function 
valued) Markov process 

y = (y t -,t>o), y t = (Y x ,x>o). 

The idea behind the approximation scheme is to first do a spatial ap- 
proximation 

V h (t)^VZ(t) = Y J ^\ 

TV 

where 7r = {r]i} i is a chosen partition of K converging to identity, and q 
are coefficients depending only on h. This approximation is uniformly 
good with respect to time, that is sup 4<T \ Vh(t) — V£{t)\ goes to zero 
when 7r grows. 
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The second step is the time dicretization, which can be performed in 
the following way : 

(4) Y m (t + At) = Y m (t) + AB t - r)iY m (t) . 

This second step is a classical one : see e.g. the book of Bouleau and 
Lepingle |J, chapter V, section B, pages 271-278, where it is stated 
that to get a precision of e, one has to consider at most a constant 
C times iV = log(^) time points. Nevertheless, when one needs to 
do the whole family of approximations (Q), one needs to control pre- 
cisely the dependence of the constant C = C(^,T) with respect to 
the space and time variables. The classical method, combined with 
Gronwall's lemma, gives an exponential upper bound for the L p -norm 
of the approximation error C(j],T) < C\(T ', p)rj exp(c2(T ', p)f] p ) . There- 
fore, we have to use another time discretization scheme for which we 
have polynomial constants (see e.g. or ||). 

The first step amounts to approximate a linear functional of an infinite 
Markov process by a linear functional of a finite dimensional Markov 
process. The secret of the success of this method is that for a large 
class of functions h, this dimension can be kept low : for a precison of 
e, the cardinal of the partition 7r(e) is roughly 

AT(e) = |7r( e )| ~c^/i log(i). 

We do not say that this approximation method is the best, but we want 
to stress the fact that the Markov property of the underlying process y 
allows an algorithm that is very economic in memory allocation. You 
only have to store N(e) values at each time step, then simulate one 
Gaussian random variable, and then update these values for the next 
step (according to rule @). 

The example of fractional Brownian motion. We can take ad- 
vantage not only of the Gaussian nature of the processes Y x , but also of 
their semimartingale nature (via Ito's formula). We have been able to 
give exact L 2 and almost sure convergence rates for fractional Brownian 
motion of index H — a — ^ e (0, |). 
Indeed, to obtain the precision 

sup\\V h (t)-V h *(t)\\ L2 <e, 

t<T 



we may choose the compact 

K(e) 

and a geometric partition 



i 



1 , 



d l-a ; ^ — J 2a+3 



n (e) = I r - m ,r- m+ \ 
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of ratio r = r(e) = 1 + c(a)y/e. 
Furthermore, if e n decreases fast enough to 0, 

(5) e n = o((togn)"$^) 

then, we have the uniform almost sure convergence 



sup 

t<T 



V h (t)-V^\t) 



almost surely . 



Theoretical results (see [0] or |§) show then that to obtain a precision 
of e = 10~ 3 , we need to consider a few hundred points. Practical re- 
sults show that forty points are enough. Furthermore, we have been 
able to check that (almost surely) the numerical approximations stick 
exponentially fast to the actual random process (and this happens uni- 
formly on compact sets), when e n goes faster to than indicated in (|5|). 

Problem II : an ergodic theorem. It is well known (Birkhoff's 
ergodic theorem) that for a stationary process U, we have the conver- 
gence, almost sure and in L 1 , 



(6) - / 0(f/( S ))d S r ^E[0(f/(O))], 

t.l Q t^oo 



as soon as <fi is an integrable Borel function: E [|0([/(O))|] < +oo. 
Of course, we first have to check the ergodicity, that is that the invariant 
cr-field is almost surely trivial, which is not so easy to do in general (see 
e.g. Rozanov |H| , page 111). 



But what we actually want to do is to answer the following question. 
Given (V(t) ; t > ) a real Gaussian process, a priori non stationary, 
such that V r = (V(t + r) ; t > ) converges in distribution to the sta- 
tionary process U. We want to show that, as soon as <ft is an integrable 
Borel function: E [|0(C/(O))|] < +oo, we have the convergence, almost 
sure and in L 1 , 



1 



(7) - / ^V(s))d sr ^E[<p(U(0))]. 
t J Q *~»-°° 

To be more precise, given h e L 2 (R + ) and (B t ,t e R) a standard 
Brownian motion, we let 

(8) V(t) = j h(t- s) dB s , and U(t) = [ h(t - s) dB s . 
We shall show that (0) holds if h is the Laplace transform 

(9) h{u) = [ e~ ux fx(dx), 



of a positive measure /j,, finite on compact sets and such that 
(10) Ml A — fi(dx) < +oo . 



x 
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2. The Markov and ergodic properties of the underlying 
infinite dimensional process 

Let B = (B t ; t > ) be a standard Brownian motion. Given x > 0, 
Y x = (Y x ; t > ) is the Ornstein-Uhlenbeck process driven by B of 
parameter x (starting from 0): 

yx = I e -x(t-s) dB _ 





It is also the solution of the stochastic differential equation 

dY x = dB t - xY t x dt . 

For every t > 0, we consider y t = (Y x , x > 0); it is a centered Gaussian 
process of covariance 

1 _ e -t(x+y) 

Fy t (x,y) = ■ (x,y>0). 

x + y 

The natural candidate for a limiting process is 3^o the centered Gauss- 
ian process of covariance 

r y 00 ( x ,y) = ^— (x,y>0). 
x + y 

Proposition 1. The process y = (yt]t> 0) is an infinite dimen- 
sional Markov process. 

Proof. The strong Markov property is quite easy to prove. Indeed, let T 
be an almost surely finite stopping time (with respect to (jF t ; t > ) the 
completed, hence right-continuous, filtration generated by the Brown- 
ian motion B). Then, B = (B t = B t +T — Bx'jt > 0) is a standard 
Brownian motion independent of (B u ,u < T), and (Y x +Tl ;t > 0) is 
the solution of the stochastic differential equation 



dY x +T = dB t — x Y x +T dt. 



Therefore, 



l Y x + [ e~ x{t - s) dB s . 
Jo 

Hence, for any bounded measurable function ip, 



I t+T — e 



E[^(Y t+T ;t >0) \F T }=Ey T [^(Y t ;t>0)], 

with 

^ v Wt ; t > )] = E [ip((Y x ' y ( x \x > 0) ; t > )] , 

and (Y t x,y ; t > ) is the Ornstein-Uhlenbeck process driven by B, of 
parameter x, starting from y. □ 
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Proposition 2. Suppose that \x is a sigma-finite measure on (0, +oo) 
such that, for a p > 1 , 



/^) SU p(*-"V-i)<+oc. 



Then, 

1. For every y e L l ([i), under F y , we have the convergence in distri- 



bution 



t^oo 



2. The only bounded invariant functions are constant. 

Proof. (1) We shall first prove the convergence of characteristic func- 
tions, since this will characterize the distribution of a possible limit 
(see ||), and then prove the tightness of the family of distributions of 
(y t , ;*>0). 

The characteristic function is defined by 

<Py t { X )=^y[e l{X ' yt) ] (Xer). 

Observe that (x, yt) is, under F y , a Gaussian random variable of mean 

/x(x)y(x)e~ xt fx(dx) i — > by dominated convergence 
and of variance 

\ 2 

x(x)Y t x jj,{dx 



E 






X (x)x(y)Kdx)fi(dy)E[Y x Y t y ] 

I _ e -t(x+y) 

X(x)x(y)ti(dx)fi(dy) ■ 

x + y 

/i(dx)fi(dy) 



. . x(x)x(y)- 

t^oo J J x + y 

This last integral is finite, Indeed, since ah < |(a 1 + 6 2 ),we can majorize 
it by 



£M -' i-2^Qr = 2 u HJ^) <+0 °' 





Therefore, (x, yt) converges in distribution to a centered Gaussian vari- 
able with variance 

H(dx)n(dy) 

x{x)x{y) : > 

x + y 

and this is easily seen to be the variance of (x, 3^oo)- 
We shall now prove the tightness of the family of laws of the processes 
(yt, t > 1). We need to introduce a little topology: we consider that for 
every t > 0, y t is a random variable taking its values in endowed 
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with the topology ^(L 1 ,^ 00 ); this is the coarsest topology for which 
the functions 

are continuous. Indeed, 



E 



c / n(dx) 



1-e 



-2xt 



2x 



< +00 



Let 4>{x) = 1 A ^. Given M > 0, Lemma ^ implies that, if 1 
1/p + 1/q, then 

\p 



is a relatively compact set. 

Let us notice that, by construction, under P y the law of y t is the law 
, under P of 3^ + ( x e~ xt y(x)). Since this last function converges 
weakly to the null function, we can suppose y = 0, without any loss in 
generality. We have the sequence of inequalities, for t > 1, 



< —E 

~ M 

< — E 

- M 



< — E 

- M 



\ 1 t I 



l-e 



2xt \ P/ 2 



2x 



P/2 



^ 

2x ) <p(x)p/i 



(f)(x)p/i 
fi(dx) 



fi(dx) 



This last quantity goes uniformly in t > 1 to zero as M — > 00 , since 
the integral is finite. 

(2) (i) First of all, we have to prove that if / is a bounded function 
depending on a finite number of coordinate, then Ptfiy) 1 — > E [/(3^oo)]- 

Indeed, we have, for a bounded measurable function g : R" — > M, and 
&l<«<nCl°°(|i) 



f(y) = g{(Xi,y), i<i<n). 



Then 



P t f(y) = E y [g((Xi,yt),l <i <n)] 

= E [^((xi.yCOe-*-) + (Xi,y t ), 1 < i < n)] . 

But we know that if X = (Xi, . . . ,X n ) is a non degenerate centered 
Gaussian vector, then we have the absolute continuity relationship 

E[g(z + X)]=E[A(z)g(X)} (z E M. n ) , 



FRACTIONAL BROWNIAN MOTION AND THE MARKOV PROPERTY 9 

with A a continuous bounded function. Applying this to ({xi, 3^t), 1 < 
i < n), for t > 0, yields 

P t /(y) = E [A((xi, 2/(.)e-'->, 1 < i < n)g(( Xl , M), 1 < i < n)] . 
Now if y p — > y weakly in then, for each i, 

(XhVpQe' 1 -) -> (Xi,2/(.)e~*-) 
and Ptf(y p ) — > Ptf(y)- Therefore, if t > 0, then P t / is a continuous 
bounded function. Since, under P„, y t i — > 34e» the law of 3^o is an 
invariant measure for the semi-group (P t ; t > 0), and we have 
lim P u / = lim P S (PJ) = E [P/(^oo)] = E [/Q^)] . 

(ii) Now, Let / : — > R be a bounded invariant function; that is, 

f(y) = P t f(y) = E y [f(y t )) {y eL\t>0). 

Observe that the sigma field T = o~(x, X £ L°°(fj,)), is generated by the 
algebra of sets depending on a finite number of coordinate 

A = A{ x ,xeL™( f i)} 

that is the smallest algebra containings the sets {y : (%, y) G B} , x £ 
P borel set of R. Therefore, the measure Py^, the law of 3^, 
being finite, there exists bounded functions (f n ,9n, n £ N) depending 
only on a finite number of coordinate such that 

fn<g<9n Mfn-9n\{yoo)} > 0. 

n— >oo 

We take the limit as £ — ► +00 in 

P t fn(y) < P t f(y) = f(y) < P t g n (y) 

to obtain 

®\fn(yoo)]<f(v)<®\9n(yoo)] 

and then let n go to 00 to conclude that f(y) = E [(3^o)] is a constant 
function. 

□ 

Lemma 3. Let n be a sigma-finite measure, and <p be a strictly positive 
4> > 0, integrable function. Then, for every M > 0, p, q > 1 such taht 
1 = l/p+l/q, 

K={f : f\I^.dfJi<M 



(fip/q 

is a relatively compact set of endowed with the weak topology 

a(L\L°°). 

Proof. Following Theorem IV. 8. 9 of Dunford-Schwarz |3| we have to 
prove that 
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1. K is a bounded subset of L 1 . 

2. If E n is a sequence of measurable sets that decreases to 0, then 



sup feK 



LJ d v 



o. 



By Holder's inequality, for / in K 



I 



^E n \f\dfi< 



f 



< M 1/p 



Lv(p) 



w 1/9 L (M) 

\ 1/5 



Taking E = f2 the whole space, we see that K is bounded in L 1 . The 
dominated convergence theorem gives that f l^ n 0d/x 

convergence is uniform for / in K. 



i — ► 0, so the 

□ 



Proposition 4. Suppose that \i is a sigma-finite measure on (0, +oo) 
such that for a p > 1 , 

J yi,(dx) sup(a;~ p/2 ,:r~ 1/2 ) < +oo . 
Then, for any measurable ip : L 1 ^) — > M such that E [IV'CVoo)!] < +°° 

7 / ^(y s )ds^E[^(y 00 )], 

t J Q t^oo 

where the convergence is almost sure and in L 1 . In particular, h(u) = 

f^°e~ ux n(du) is in L 2 (0,+oo); hence, if we let N denote a standard 

i 

Gaussian random variable, a = ||^|| L 2( M+ ) = (f^°h(u) 2 du) 2 , and 

V h (t) = [ h{t-u)dB u , 
Jo 

then, for any measurable ip such that E [^(iV) |] < +oo ; for any y e 
under F y , 



1 f* 1 

-V h (s)) ds i — > 

a t^oo 



t 



E[0(iV)] 



o 



where the convergence takes place almost surely and in L 1 . 

Proof. The second part of the proposition is an obvious consequence of 
the first part: we only need to take ip(y) = (f)(f y(x) d/i(x)) and to use 
Fubini-Tonelli's theorem 



POO POO 

/ h(u) 2 du = du 
Jo Jo 




d/i(x) dji{y) 
d/i(x)dfj,(y) 



-yu 



d/x(y) 



x + y 
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To prove the first part, we consider the invariant sigma- field X, 

F = jz : t £ R+ -»• zt £ L\fi), ijf ^(z s ) ds h— > E [^(^oo) I X] J , 

and /(?/) = Pj/((t — » 3^t) ^ -F). Birkhoff's ergodic theorem implies that 
E[/(D4o)] = (and that the convergence also takes place in L 1 ). 
The invariance of / is an easy consequence of the Markov property 
of (y t ;t > 0) since for all t > 0: F — 9^{F) where Q t denotes the 
shift operator. Proposition H implies that / is a constant, and thus, for 
all y e L 1 ^), f(y) = 0. 

All we need to show now is that the invariant sigma-field X is almost 
surely trivial. Let Z be a bounded X measurable random variable. We 
have, for all t, Z = Z o 9 t and therefore the function f(y) = ~E y [Z] is 
constant since 



(11) E[Z\F t )=E[Zoe t \F t ] = f(y t ). 

Taking expectations under F y yields 

P t f(y) = E y [f(y t )] = f(y). 

Proposition |^ implies that / is constant, say f(y) = c for every y £ 
L 1 (//). Injecting this in relation ( |TT| ) yields, that under F y 

E[Z\J 7 t ] = c (t > 0) . 

Taking limits as t — > +00 gives finally that Z = c a.s. 

□ 
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